FINITE ELEMENT EXTERIOR CALCULUS FOR PARABOLIC 

PROBLEMS 

DOUGLAS N. ARNOLD AND HONGTAO CHEN 

Abstract. In this paper, we consider the extension of the finite element exterior calculus 
from elliptic problems, in which the Hodge Laplacian is an appropriate model problem, to 
parabolic problems, for which we take the Hodge heat equation as our model problem. The 
numerical method we study is a Galerkin method based on a mixed variational formulation 
and using as subspaces the same spaces of finite element differential forms which are used for 
elliptic problems. We analyze both the semidiscrete and a fully-discrete numerical scheme. 



1. Introduction 

In this paper we consider the numerical solution of the Hodge heat equation, the parabolic 
equation associated to the Hodge Laplacian. The initial-boundary value problem we study 



is 

(1.1) ut + {dS + Sd)u = f mnx{0,T], 

(1.2) tT{i.u) = 0, tT{i^du) = onOnx (0,T], 

(1.3) m( ■ , 0) = Mo in Q. 



(We could consider other boundary conditions as well.) Here the domain Q C M" has a 
piecewise smooth, Lipschitz boundary, the unknown m is a time dependent differential Re- 
form on Q, Ut denotes its partial derivative with respect to time, and d, S, -k, and tr denote 
the exterior derivative, coderivative, Hodge star, and trace operators, respectively. 

The numerical methods we consider are mixed finite element methods. These are based 
on the mixed weak formulation: find (a, m) : [0, T] — )■ HA.^~^ x Hh!', such that 'u(O) = uq 
and 

(1.4) (a,r) - {dT,u) = 0, r G HK^'^, t G (0,T], 

(1.5) {ut,v) + {d(7,v) + {du,dv) = {f,v), v e //A^ t e (0,T], 

(The notations are explained in the following section.) Notice that, unlike in the elliptic 
case, the harmonic forms need not be explicitly accounted for in the weak formulation. 
The well-posedness of the mixed formulation (11. 4p . (11. 5p is established in a precise sense in 
Theorem 14.41 

In the simplest case of 0-forms {k = 0), the differential equation (II. ip is simply the 
heat equation, Ut — Au = /, and the boundary condition (II. 2p is the Neumann boundary 
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condition, du/dn = 0. Moreover, in this case the space HA'' ^ vanishes, and the weak 
formulation f ll.4p - fll.5p is the usual (unmixed) one: u : [0,T] — )• H^{Q) satisfies 

{ut,v) + (gradM,gradt;) = {f,v), v G H\n), t G (0,T]. 

In this case, the numerical methods and convergence results obtained in this paper reduce 
to ones long known [HI [10] . In the case of n- forms, the differential equation is again the heat 
equation, although the natural boundary condition is now the Dirichlet condition u = 0. In 
the case of ra-forms, the weak formulation seeks cr g iJ(div), u E Li^ such that 

((T,r) - (divr,n) = 0, r G if(div), {ut,v) - (diva,?;) = (/,i;), v G te (0,T]. 

This mixed method for the heat equation was studied by Johnson and Thomee in [S] in two 
dimensions. Recently, Hoist and Gillette [S] have studied this mixed method in n-dimensions 
using a finite element exterior calculus framework (in their work they consider hyperbolic 
problems as well). 

For A; = lor2inn = 3 dimensions, the differential equation (11. ip is the vectorial heat 
equation, 

Ut + curl curl u — grad div u = f. 

The weak formulations (ll.4p -( IT3]) for k = 1 and 2 correspond to two different mixed formu- 
lations of this equation, the former using the scalar field a = div u as the second unknown, 
the latter using the vector field a = curlw. For k = 1, the boundary conditions (II. 2p . which 
are natural in the mixed formulation, become u - n = 0, curlw x n = 0, while for k = 2 these 
natural boundary conditions are u x n = 0, divu = 0. This vectorial heat equation arises, 
for example, in the linearization of the Ginzburg-Landau equations for superconductivity 
[7], and is related to the dynamical equations of Stokes and Navier-Stokes fiow (see, e.g.. 

To discretize (II. 4p . (II. Sp . we utilize the two main families of finite element differential 
forms, the P^-A^ and P~A'^ spaces. Between them they include lots of the best known 
families of finite elements on simplicial meshes [2l Section 5]. We give both semidiscrete 
and fully discrete schemes, and the corresponding convergence analysis. Convergence rates 
under different norms are shown in our final results (see Theorem l5.3l and Theorem l6.3l below). 
These achieve the optimal rates allowed by the finite element spaces provided some regularity 
assumptions are satisfied. These results also reveal the relation between convergence rates 
under different norms and the regularity of the exact solution. 

The outline of the remainder of the paper is as follows. In Section [21 we review of basic no- 
tations from finite element exterior calculus, including the two main families of finite element 
differential forms, the P~A'^ and P^A'^ families, and some of their properties. In Section [3l 
we apply the elliptic theory to define an elliptic projection which will be crucial to the error 
analysis of the time-dependent problem, and to obtain error estimates for it. In Section [H 
we turn to the Hodge heat equation at the continuous level and establish well-posedness 
of the mixed formulation. We then give a convergence analysis for the semidiscrete and 
fully discrete schemes in Sections \5\ and [Q respectively. Finally, we present some numerical 
examples confirming the results. 
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2. Preliminaries 

We briefly review here some basic notions of finite element exterior calculus for the Hodge 
Laplacian. Details can be found in [11 § 2] and [21 §§ 3-4] and in numerous references given 
there. 

For Q a domain in and k an integer, let L^A'^ = L^A'^(f2) denote the Hilbert space of 
differential /c-forms on Q with coefficients in L^. This is the space of functions on Q with 
values in Alt'^ M", a finite dimensional Hilbert space of dimension (understood to be if 
k < or k > n). The Hodge star operator -k is an isometry of Alt'^ M"" and Alt"~''"]R", and 
so induces an isometry of L'^A^ onto L^A""*^. The inner product in L^A'^ may be written 
{u, v) = uA-kv, with the corresponding norm denoted We view the exterior derivative 
d = d^ as a unbounded operator from L^A'^ to L^A'^"'"^. Its domain, which we denote iJA^(r2), 
consists of forms u G L^A'^ for which the distributional exterior derivative du belongs to 
^2^fc+i_ Assuming, as we shall, that Q has Lipschitz boundary, the trace operator tr = tig^ 
maps HA''{il) boundedly into an appropriate Sobolev space on dil (namely H~^^'^A''{dD.)). 
The coderzvatzve 6 is defined as (-l)"('^+^)+i*ci^ : H*A'' -> H*A''-\ where H*A'' := *i7A"-^ 
The adjoint d* = d^ of d^~^ is the unbounded operator L'^A'' — )■ LF'A^'^ given by restricting 
5 the domain of d* , 

D{d*) = H*A^ ■.= {ue H*A^ I tr *M = }. 

We denote by 3^^ and 3^ the null spaces of d^ and ci^, respectively. Their orthogonal com- 
plements in LF'A^ are QS^ and fB'^, the ranges of and d^~^, respectively. The orthogonal 
complement of inside 3^ is the space of harmonic forms 

i^'^ = 3^= n 3fc = { G HA^iSi) n h*a\q) \du = o, d*u = o }. 

The dimension of i^'^ is equal to the kth Betti number of Q, so i^i'^ = for /c 7^ if f2 is 
contractible. The Hodge decomposition of L'^A^ and of HA'' follow immediately: 

(2.1) L^A'' = ^'' 

(2.2) HA'' = ^'' ®S^'' 

where 3'^"'" = HA'' fl denotes the orthogonal complement of 3'^ in HA''. 

The Hodge Laplacian is the unbounded operator L = dd* + d*d : D{L) C L^A'^ — )■ L'^A'' 
with 

D{L) = {ve HA'' n H*A'' I d*v G HA''-\ dv G H*A''+^ }. 

The null space of L consists precisely of the harmonic forms S)''. 

For any / G L'^A'', there exists a unique solution u = Kf G D{L) satisfying 

Lu = f (modii), u ± S^'', 

(see [21 Theorem 3.1]). The solution u satisfies the Hodge Laplacian boundary value problem 

{d5 + 5d)u = f — P^f in fi, trku = 0, tr^trrfu = on d^l, 

together with side condition u S^'' required for uniqueness. The solution operator K is a. 
compact operator L'^A'' — )■ HA'' fl H*A'' and a fortiori, is compact as an operator from L^A^ 
to itself. 
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Now we consider the mixed finite element discretization of the Hodge Laplacian boundary 
value problem, following [2]. This is based on the mixed weak formulation, which seeks 
a e HA^~^, u G HK^, and p e S)'' such that 

(cr, r) - (dr, n) = 0, r G i/A^~\ 

(da, v) + {du, dv) + (p, v) = {f,v), v e HA'', 

{u,q) = 0, qe^\ 

It admits a unique solution given hy u = Kf, a = d*u, p = Pf,f. We discretize the mixed 
formulation using Galerkin's method. For this, let A^~^ and A^ be finite dimension subspaces 
of HA.^~^ and HK^, respectively, satisfying dA^^ C A^. We define the space of discrete 
harmonic forms S)i\ as the orthogonal complement of := dA^"^ inside 3^ := 3 H A^. This 
immediately gives the discrete Hodge decomposition 

where 3^"*" is the orthogonal complement of 3/^ inside A|. 

The Galerkin method seeks ah E A^~^, Uh G A^, ph G Sj^ such that 

{(rh,T) - {dT,Uh) = 0, rGA^-\ 
(2.3) {dah, v) + {duh, dv) + (p/,, t;) = {f,v), v e A^, 

(u,„g) = 0, gGfit 

For the analysis of this discretization, we require the existence of a third space A^^^ C 

HA''+^ which contains dA^, so that A^"^ 4 A^ 4 A^+^ is a subcomplex of the segment 

HA'' A HA^^^ of the de Rham complex. Further we require that there exists 
a hounded cochain projection, i.e., bounded linear projection maps vr^^ : HA^ — A;^, j = 
k — 1, k, k + 1, such that the diagram 



(2.4) 



.fc_l d .fc , Afc+1 



commutes. A key result of the finite element exterior calculus is that, under these assump- 
tions, the Galerkin equations (12. 3p admit a unique solution and provide a stable discretiza- 
tion. 

Another important aspect of the finite element exterior calculus is the construction of 
finite element spaces A^ which satisfy these hypothesis, i.e., which combine to form de Rham 
subcomplexes with bounded cochain projections. Let there be given a shape regular family 
of meshes Th with mesh size h tending to 0. For each r > 1, we define two finite element 
subspaces of HA'', denoted VrA'^iTh) and V~A''{Th)- For k = 0, these two spaces coincides 
and equal the degree r Lagrange finite element subspace of H^{Q). For k = n, V~A'^{Th) 
coincides with Vr-iA"'{Th), which may be viewed as the space of all piecewise polynomials of 
degree at most r — 1, without inter-element continuity constraints. However, for < k < n, 

Vr-iA^iTh) C V^A^iTh) C VrA^iTh). 
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For stable mixed finite elements for the Hodge Laplacian, we have four possibilities (which 
reduce to just one for k = and to two for k = 1 or n): 

or 

Vr^iA'^iTh) (ifr>i: 

As the auxihary space, if = V^A^iTh), we take A^+^ = V~ A''+^ {%) , while if = 
Vr-iA'^iTh), we take A^+i = V^^^A'^+^Th). 

For this choice of spaces, it is known ([H § 5.4], [21 § 5.5], [5]) that there exist cochain 
projections as in (12. 4p for which vr^ : L'^A^ — t- A^ is bounded in L^A-' uniformly with respect 
to h. In particular, this implies that there is a constant C independent of h such that 







(VrA^-HTh)) 




(2.5) 












[v-A^~\%)\ 





(2.6) \\u-'Kiu\\<C mf \\u-vl ueL^A^. 
Moreover, we have the approximation estimates 

(2.7) \\u-7iiu\\ <Ch'\\u\\s, 0<s< 
Note that we use ||m|L as a notation for the Sobolev norm ||M||/i-sAj 




Ai = p-A^(r.), 

1, Ai = VrA^{Th). 



3. Elliptic projection of the exact solution 

As usual, we shall obtain error estimates for the finite element approximation to the 
evolution equation by comparing it to an appropriate elliptic projection of the exact solution 
into the finite element space. In this section we define the elliptic projection and establish 
error estimates for it. 

Given any u E D[L), the elliptic projection of u is defined as [ah, Uh-,Ph) G A^~^ x x i^^, 
such that 

(3.1) {ah,T)-{dT,Uh) = Q, teAI-\ 

(3.2) {ddh, v) + {duh, dv) + {ph, v) = {Lu, v), v e A^, 

(3.3) {Uh,q) = {u,q), q E Sj^ 

By Theorem 3.8 of [2] there exists a unique solution to f l3.ip -f lX^ . Now we follow the 
approach of [2] to derive error estimates. To this end, we introduce some notations. First, 
let Ps)^ : L^A^ i^^ denote the L^-projection. From (13. 3p . Psj^Uh = -Pfih^. Moreover, from 
[21 Section 3.4], 

Ph = P^^n^Lu) = Pfijrfa), 

where a = d*u, the last equality holding because d*du G 53^ ± 3^, but i^^ C 3^ C 3^^- 
Next, define f3 = (3^, fi = jjf^, and rj = rj^hy 

/3 = ||(/ - 7r/jK||£(i2Afc,L2Afc), /i = ||(/ - 7r/j)P^fe||£(2,2Afc,L2Afc), 

?7 = maxmax[||(/ - TXh)dK\\nL2p^k~j L2f^k-j+i), - Tih)d* K\\ciL2[^k+j L2[^k+j-i)]. 
i=o,i 
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From (12. 6p and the compactness of : L^A'^ — )■ Hh^ fl H*A'', we conclude that 77, /3, /i — )■ 
as /i — i- 0. Assuming regularity for the Hodge Laplacian (by which we mean both that 
\\Kfh < Cll/llo for all / G L^A'^ and that Sj'' C H^), then we have 

(3.4) 7] = 0{h), P,fi = 0(/l"^'°(2,r)) 

for any of the choices of spaces in (12.51) . Note that /3 = 0{h?) except in the case r = 1 and 
so A| = A^ 

Finally, we denote the best approximation error in the norm by 
E{w) = inf \\w - v\\, w e L^A'', k = 0,. . .,n. 

We are now ready to give the error estimates for the elliptic projection. 

Theorem 3.1. Let u G D{L) and let {ah,Uh) be defined by (I3.ip - (l3.3p . Then we have 



(3.5) \\d{a-&,)\\<CE{da), 

(3.6) ||a - ahW < C{E{a) + r]E{da)), 

(3.7) WphW < CfxE{da), 

(3.8) \\d{u-Uh)\\ <C{E{du)+riE{da)), 



(3.9) \\u -Uh\\< C{E{u) + E{P^u) + 7][E{du) + E{a)] + (r/^ + p)E{d(T) + fxE{P<Bu)). 

Proof. This is essentially proven in [2], except that there it is assumed that u -L Sj and 
ith J- Sjh- To account for this difference, let m = m — Pf,u and Uh = Uh — Pfih^h- Then (13.11) 
and (13.21) continue to hold with u and Uh replaced by u and m^, respectively, and, in place 
of (13. 3p . we have 

{uh, g) = 0, g G fit 

Application of Theorem 3.11 of [2] (with f = Lu and p = 0) then gives the (I3.5ti3.8p . and, 
instead of (13. 9p . we get 

\\u - Uh\\ < C{E{u) + 7][E{du) + E{a)] + {rf + /3)E(rfa) + iJiE{P^u)) 

< C{E{u) + ^(Py,M) + r][E{du) + E((t)] + (r/2 + l3)E{da) + fiE{P^u)). 

Thus ||tt — tt/j|| is bounded by the right-hand side of (13.90 . and, to complete the proof, it 
suffices bound Pf^u — Psj,^Uh by same quantity. Now 

Ps^U - Ps^i^Uh = Pf,U - Pfj^U = (/ - PfjJPs^U - Ps)i^{u - PfjU), 

For the first term on the right-hand side, we use [21 Theorem 3.5] and (12. 6 p to get 

- P^jP,,u\\ < \\iI-7rh)P^u\\<CE{P^u). 

To estimate the second term, we use the Hodge decomposition (12. 2p to write u—Ps^u = Ub+u± 
with Ub G ^^,u± G S''"^- Since C 3'', Pfjh'^^ = 0; since iVhUb G ®t Psj^'^hUb = 0. 
Hence P^^{u — Pfju) = PfjiX^ — i^hjUb- We normahze this quantity by setting 

q = P^,Su - P^u)/\\P^,Xn - P^u)\\ G fit 

Then P^-^q G i^, and, by [2], Theorem 3.5], ||g — P^^q\\ < — 7r/i)Pf,g|| < /i. Therefore, 

\\Pf)hk^ - P^-)U)\\ = {PsjiSu - Ps)u),q) = {PsjiSi - Trh)ub,q) = ((/ - 7rh)ub,q). 
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Now (J — HhjUb G ^B^, and so is orthogonal to S). Thus 

((/ - Tih)ub,q) = ((/ - 7rh)ub,q- Ps^q) < - Tih)ub\\\\q - i^gll < CiiE{P<su), 
by (12. 6p . Combining these results, we get 

\\P,,u - P,,^Uh\\ < C[E{P,,u) + ^lE{P^u)l 
completing the proof of the theorem. □ 

Assuming sufficient regularity of u and a = d*u, we can combine the estimates of the 
theorem with the approximation results of (12. 7p to obtain rates of convergence for the elliptic 
projection. The precise powers of h and Sobolev norms that arise depend on the particular 
choice of spaces in f E3]) . For example, if we take A^"^ = VrA''~^{Th), then we can show 
the optimal estimate \\a — ah\\ < C/;.'''"'"^||cr||r-+i, but, if A^~^ = V~ h^~^{Th)-i then clearly we 
can only have ||cr — cr/i|| = 0{h'^). Rather than give a complicated statement of the results, 
covering all the possible cases, in the following theorem and below we restrict to a particular 
choice of spaces from among the possibilities in (12. 5p . Moreover, we assume r > 1, since 
the case r = 1 is slightly different. However, very similar results can be obtained for any of 
the choices of spaces permitted in (12.50 . including for r = 1, in the same way. Finally, we 
introduce the space 

= {ueH'\ P^u e H', P^u G H'-^ }, 

with the associated norm 

since it will arise frequently below. 

Theorem 3.2. Assume regularity for the Hodge Laplacian, so (13.41) holds and suppose 
that we use the finite element spaces A^~^ = V~ A^~^{Th) and A^ = V~A''{Th) (so that 
the auxilliary space is A^"*"^ = V~A^^^{Th)), for some r > 1. Then we have the following 
convergence rates for the elliptic projection: 

\\d{a-ah)\\ < Ch''\\da\\r, 
\\a — ay^W < C^'WaWr, 
\\ph\\ < Ch''\\da\\r^2-, 
\\d{u - Uh)\\ < Ch'{\\du\\r + 

\\u — Uh\\ < CK'WuWht- 

4. Well-posedness of the parabolic problem 

We now turn to the Hodge heat equation. In this section we demonstrate well-posedness 
of the initial-boundary value problem (II. 4p . (II. 5p . The key tool is the Hille-Yosida-Phillips 
theory as presented, for example, in[3] and 

We begin by showing that the Hodge Laplacian is maximal monotone, or, equivalently, in 
the terminology of [1], that its negative is m-dissipative). This is the key hypothesis needed 
to apply the Hille-Yosida-Phillips theory to the problem (ll.ip - (ll.3p . 
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Theorem 4.1. The Hodge Laplacian L is maximal monotone. That is, it satisfies 

{Lv,v)>0, \/veD{L), 
and, for any f G L'^A^, there exists u G D[L) such that u + Lu = f . 

Proof. For any v G D{L), {Lv,v) = {dv,dv) + {d*v,d*v), so the monotonicity inequahty 
is obvious. Now, for any / G L^A'', the Riesz representation theorem furnishes a unique 
u G HA'' n H*A'' such that 

(4.1) {du,dv) + {d*u,d*v) + {u,v) = {f,v), veHA''r]H*A\ 

We shall show that this u belongs to D{L), from which it follows immediately that u+Lu = f. 

To show that u G D{L), we must show that du G H*A'''^^ and d*u G ifA^~^. From 
(14. ip . f — u is orthogonal to S^'', so, using the Hodge decomposition of L^A'^, we may write 
f -u = df\ + d*f2 with /i G HA''-^ n and /a G if^A'^+i n Then 

(/ -u,v) = {df\ + d*f2, v) = (A, d*v) + (/2, A^), G i/A'^ n if*A^ 
Combining with fl4.ip . we get 

(4.2) {du - /2, rft;) + {d*u - /i, = 0, ve HA'' n ^*A^ 

Now du, f2 G *B''+\ so there exists v G 3'''^ = HA'' fl 23^ such that dv = du - f2. Choosing 
this V in (14.20 . we find du = f2 ^ H*A''~^^, as desired. Similarly d*u = /i G HA''"^. □ 

Since L is maximal monotone and self-adjoint, we obtain the following existence theorem. 
This is proved in [3] in Theorems 3.1.1 and 3.2.1 for / = and Uq G L^A*^, and in Proposi- 
tion 4.1.6 for general / and uq in D{L). Combining the two results by superposition, gives 
the theorem. 

Theorem 4.2. Suppose thatuo G L^A'^ and f G C([0, T]; L^A'^) are given and that either f G 
L\{0,T);D{L)) or f e W^^\{0,T); L^A'') . Then there exists a unique ue C([0, T]; L^A'^) n 
C((0,T];D(L)) nCH(0,T];L2A'=), such that 

Ut + Lu = f on Q X (0, T], u{0) = Uq. 

If further, uq e D{L), then u G C{[0,T]; D{L)) n C\[0,T]; L'^A''). 

We denote by S{t) : L^A^' — L'^A'' the solution operator for the homogeneous problem 
(/ = 0), so u{t) = S{t)uo solves Ut + Lu = 0, m(0) = uq. Then S{t) is a contraction in L'^A'' 
for aU t G [0,T], i.e., ||^(t)|| < 1, and S{t) commutes with L on D{L) (Theorem 3.1.1 of j3|). 

We can measure the regularity of the solution (for general /) by using the iterated domains 
defined by D{V) = {u e D{L^~^) \ V^^u G D{L) }, / > 2. The next theorem show that if / 
is more regular, then the solution is also more regular. 

Theorem 4.3. Suppose that in addition to the hypotheses of Theorem \4.^ we have that f 
belongs to C((0, T]; D(L)) n ^^((0, T); ^(L^)). Then 

(4.3) ueC\{0,T];D{L)). 
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Proof. If / = 0, then [3l Theorem 7.7] imphes that 

ueC\{0,T]-D{L^)), 

for all k,l > 0. Therefore, it is sufficient to treat the case Uq = 0, which we do using 
Duhamel's principle. By Proposition 4.1.6 of the solution is given by 

u{t) = [ S{t- s)f{s)ds 
Jo 



in this case, and, assuming that / satisfies the hypotheses of Theorem 14.21 

u e C([0, T];D{L)) n C\[0,T]- L^K^). 
Now / G Li((0,r);D(L2)), so 



L^u(t) = [ S{t- s)L^f{s) ds, 
Jo 



by the commutativity of S{t — s) and L. Since S{t — s) is a contraction in L^A'^, this 
implies that u G C([0, T]; Z}(L^)) and so Lu G C {[0 , T]; D (L)) . Since we also assume that 
/ G C((0,T]; D{L)), ( 14. 3p follows immediately from the equation ut = f — Lu. □ 

Next we show that the solution u guaranteed by Theorem 14.21 together with a = d*u, 
is a solution of the mixed problem (11. 4p . (11. 5p . Since u G C((0, T]; a = d*u G 

C{{0,T];HA^-^) and (O]) holds. Clearly 

{ut, v) + {Lu, v) = {f,v), ve L^K\ t G (0, T]. 

Since u G C((0, T]; £)(L)), we have 

{Lu, v) = {dd*u, v) + {d*du, v) = {da, v) + {du, dv), v G HA'', t G (0, T]. 

Combining the last two equations gives (II. 5p . 

We are now ready to state the main result for this section. 

Theorem 4.4. Suppose that uq G L'^A'' and f G C([0, T]; L^A'^) are given and that either 
f G L^{{0,T);D{L)) or f e W^'\{0,T); L'^A'') . Then there exist unique 

a G C((0,T];i7A^"i), u G C([0, T]; L^A'^) n C((0, T]; n C^((0, T]; L^A'^), 

satisfying the mixed problem (11.40 . (11.51) and the initial condition m(0) = Uq. If, moreover, 
the hypotheses of Theorem \4.3\ are satisfied, then (14. 3 p holds. 

Proof. We have already established existence. For uniqueness, we assume / = and take 
r = cr in (ll.4p and v = u in (II. 5p . to obtain 

~\\ur = -\wr-\\dur<o. 

Therefore is decreasing in time, so if n(0) = 0, then u = 0. Finally, (ll.4p then implies 
that a = 0. □ 
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5. The semidiscrete finite element method 

The semidiscrete finite element method for the Hodge heat equation is Galerkin's method 
applied to the mixed variational formulation fll.4p . (11. 5p . That is, we choose finite element 
spaces A^~^ and as in (12. 5 p for some value of r > 1, and seek {ah, Uh) E C([0, T]; A^^"^) x 
C^([0,T]; A^), such that ^^(O) = u^, a given initial value in A^, and 

(5.1) (a^,r) - (rfr,n,) = 0, r G A^\ t E (0,T], 

(5.2) {uh,t,v) + {dah,v) + {duh,dv) = {f,v), vE A^ t E (0,r]. 

In this section we shall establish convergence estimates for this scheme. 

We may interpret the semidiscrete solution in terms of two operators, c/^ : A^ — A^~^ 
and Lh : A^ — )■ A^, which are discrete analogues of d* and L, respectively. For f G A^, 
dj^v E A^~^ is defined by the equation 

{dlv,T) = {v,dT), TEAI-\ 

and the discrete Hodge Laplacian L/^ : A^ — )■ A^ is given by Lh = d\d + dd\. The following 
characterization is then a direct consequence of the definitions. 

Lemma 5.1. The pair {ah,Uh) E C([0, T]; A^^^) x CH[0,T];A^) solves and if 

and only if Uh(t) E C^([0,T];A|) solves 

(5.3) Uh,t + LhUh = Phf, 0<t< T, 
where Ph is projection of f onto A^, and ah = d\uh- 

From the theory of ordinary differential equations, there exists a unique solution Uh E 
C^([0,T]; A^) solving the ODE (15. 3p and taking a given initial value. Letting ah = d\uh, we 
obtain a unique solution to the semidiscrete finite element scheme (15. ip . (15. 2p . 

Remark. The formulation (15.31) is useful for theoretical purposes, but is typically not imple- 
mented directly, rather only implicitly via the mixed method. This is because the operator 
d*^ is not local. Even if the finite element function v is supported in just a few elements, d^v 
will generally have global support. 

Next, we turn to the convergence analysis. In Proposition 15 .21 we shall give error estimates 
for the difference between the semidiscrete finite element solution and the elliptic projec- 
tion of the exact solution of the evolution equations. Combining these estimates with the 
estimates from Section [3] for the elliptic projection gives error estimates for the semidiscrete 
finite element method, which we present in Theorem 15. 3[ 

Assume the conditions of Theorem 14.31 hold, so the exact solution 

u E C([0,T];L2a'=) nCi((0,T];D(L)). 

For each t > 0, we can then define the elliptic projection of u{t) and of Utit); see (13. ip - 
(13. 3p . Writing {ah{t),Uh{t),phit)) for the former, it is easy to see that its time-derivative, 
{o'h,t,Uh,t,Ph,t), is the elliptic projection of Ut. From Theorems 13.11 and 13.21 we obtain error 
estimates, such as 

(5.4) \\ut - Uh,t\\ < C{E{ut) + E{P^ut) + r^[E{dut) + E{at)] 

+ (r/2 + l3)E{dat) + ^iE{P^ut)) < Ch'\\ut\\Hr, 
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with the last inequahty holding for the choice of spaces made in Theorem 13.21 (and similar 
results holding for the other allowable choices of spaces). Now, from (13. ip . 

(5.5) (a,,, r) - {dr, n,,) = 0, r G t E (0, T], 
and, substituting Lu = —Ut + / into ( 13. 2p . 

(5.6) {uh,u v) + {dah, v) + {duh, dv) = {uh,t - Ut, v) + (/, v) - (ph, v). 
Define 

S/j = (Th, - CT/i, Uh = Uh-Uh, 

the difference between the elliptic projection and the finite element solution. Subtracting 
(15. ip and (15. 2p from (15. 5p and (15. 6p . respectively, gives 

(5.7) (S,,, r) - {dr, Un) = 0, r G < t < T, 

(5.8) {Uh^u v) + {dT.h, v) + {dUh, dv) = {uh,t - ut - ph, v), v e AI 0<t<T. 

We shall now use these equations to derive bounds on S/^ and Uh in terms of Uh^t — Ut and 
Ph, for which we derived bounds in Section [31 In the remainder of the paper, we adopt the 
notation || ■ 11^00(^2) for the norm in L°°{0,T; L^A^(f2)) and similarly for other norms. 

Proposition 5.2. Assume uq G D{L). Then 

\\Uh\\L°°(L2) + ||S/i||l2(2,2) + \\dUh\\L^(L'2) < C{\\Uh{0)\\ + \\Uh,t - Ut- P/i||l1(L2)), 

||Sft||Loo(2,2) + II rfS/i 11^2(2,2) < C {\\dlUh{0)\\ + \\Uh,t - Ut - Ph\\L^{L2)). 

Proof. By Theorem 1121 u G C([0, T]; D(L)) n ^^([O, T]; L^A'^). For each t G (0,T], take 
r = S/,(t) G A^-^ in (JO) and v = Uh{t) G A| in (JSlHD, and add to obtain 

(5.9) ^l"^""' + ll^'^ll' + "'^^'^11' = ^^'^'^ " " 
which implies 

^iWhW^ < 2||n/,,t - -]3?,||||f//,||. 
dt 

Taking t* G |0,T] such ||f//i||L°°(L2) = ||f//j(t*)||, and integrating this inequality from to t* 
gives 

||f/h(t*)|P < ||t^h(0)|P + 2\\Uh,t - Ut - Ph\\L^iL^)\\Uh\\L°°(L2), 

whence 

(5.10) ||f//i||L-'(L2) < \\Uh{0)\\+2\\Uh,t-Ut-ph\\mL^), 

which gives the desired bound on Uh- To get the bound on S/^ and dUh, integrate (15. 9p over 
t G [0, T]. This gives 

||Sh||i2(i2) + \\dUh\\l2^L2) < ^\\Uh{0)\\^ + \\Uh\\L°°(L2)\\uh,t - Ut -P/i||li(l2), 
and so, by (I5.10p . 

||S/jL2(i2) + \\dUh\\L^L^) < C{\\Uh{0)\\ + \\Uh,t - Ut -P/j||l1(l2)), 

which completes the proof of the first inequality. 
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To prove the second inequality, we differentiate (15. 7p in time and take r = S/^ G A^^ ^, and 
then add to (15.81) with v = d'Zh £ (here we use the subcomplex property (iA^~^ C A^). 
This gives 

- — ||Sftf + ||c/S/,f = {uh,t-Ut-ph,dT.h). 

By integrating in time, first over [0, t*] with t* G [0, T] chosen so that ||Eft||ioo(/^2) = ||E;j(t*)||, 
and then over all of [0,T], we deduce that 

II S/i 11^00(^2) + WdEhWl^L^) < C'(||Sh,(0)|| + \\Uh,t - Ut -Ph\\L^{L^))- 

Finally, we note from (15.11) and (13. ip that S/^ = d'lUh, and so complete the proof. □ 

Now suppose, for simplicity, that we choose the initial data m° to equal the elliptic pro- 
jection of Uo- Then Uh{0) = and the right-hand sides of the inequalities in Proposition 15.21 
simplify. Bounding them using Theorem 13.21 and (15.41) we get, for the choice of spaces indi- 
cated in the theorem, 

||t^/i||L°°(L2) + ||S/j||i2(i2) + \\dUh\\L^{L^) < Ch'^{\\Ut\\Li(^Hr-^ + II 11^1(^,-2)), 
||Sh||L°o(2,2) + \\dEh\\L^(L'2) < Ch^dlUtWi^^H^) + \\dd*u\\L2(^Hr-2-^). 

Combining these estimates with the estimates in Theorem 13.21 for the elliptic projection, 
we obtain the main result of the section. 

Theorem 5.3. Suppose that, in addition to the hypotheses of Theorem \ 3.2\ and \4.3[ Uq G 
D{L). Let {a,u) be the solution of (T^, U.5\) satisfying U.3\) . and {ah,Uh) the solution of 
( 15. 1]) . ( 15. i^) with the spaces selected as in Theorem \3.^ and Uh{0) chosen to be equal to the 
elliptic projection of Uq. Then, we have the following error estimates for the semidiscrete 
finite element method: 

\W - (^h\\L2{L2) < Ch'WlUtWL^Hr) + ||(i*n||i2(^r)), 
||cr — 0"/,,||l°°(L2) < Ch'^{\\Ut\\L2(^Hr) + ||(i*n||Loo(^r)), 

\\d{a - (Th) 11^2(^2) < Ch''{\\ut\\L2(Hr) + \\dd*u\\L2(Hn)^ 

\\U - M/i||l°°(L2) < C/l''(||M||i«=(5r) + \\ut\\ LI 

\\d{u - n/,) 11^2(^2) < Ch''{\\ut\\mH-) + \\du\\L2(^Hn + \\dd*u\\L2(^H-^))- 

6. The fully discrete finite element method 

If we combine the semidiscrete finite element method with a standard time-stepping scheme 
to solve the resulting system of ordinary differential equations, we obtain a fully discrete 
finite element method for the Hodge heat equation (11.41) . (11.51) . For simplicity, we use 
backward Euler's method with constant time step At = T/M. We may choose any of the 
pairs of finite element spaces indicated in (12. 5p for any value of r > 1, but, as above, for 
simplicity we restrict ourselves to the choice A^~^ = V~A''~^{Th) and A^ = V~ K^{Th) with 
r > 1, the results for the other cases being simple variants. The fully discrete method seeks 
cr^ G K^^^vll G A^, satisfying the equations 

(6.1) K,r)-(dr,<) = 0, rGA^i 



J. 
At 



(6.2) ( " ,v) + {dal,v) + {dul,dv) = {f{t-),v), v E At 
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for 1 < n < M. It is easy to see that this hnear system for is invertible at each time 

step. We initiahze by choosing u\ G A^. We also define cr° G A^~^ so that (16.11) holds for 
n = 0. 

Next, we turn to the convergence analysis. We first obtain error estimates for the difference 
between the fully discrete finite element solution and the elliptic projection of the exact 
solution of the evolution equations. These are stated in fl6.10p and (16. lip . Combining these 
estimates with the estimates from Section [3] for the elliptic projection, we obtain the error 
estimates for the fully discrete finite element method presented in Theorem 16. 3[ 

The analysis is similar to that for the semidiscrete finite element method, but with some 
extra complications arising from the time discretization. Let {a^^vJ^^p^ be the elliptic 
projection of m" = uit^). 

Now, from (13.10 

(6.3) {al, r) - {dr, ul) =0, r G < n < M, 
and, from (13.21) and the equation ut + Lu = /, 

(6.4) ( ' , v) + {dal v) + {dul dv) = ( ' - <, v) + {f\ v) - {pi, v) 

« - n") - (m;J-i - n'^-i) , 



Set 



At ' ' ^ At 

w G At l<n<M. 



spn ±.n n jjn '-n n 



the difference between the elliptic projection and the finite element solution at each time 
step. Subtracting (16. ip and (16. 2 p from (16. 3p and (16. 4p . respectively, gives 

(6.5) (S^, r) - (rfr, U^) = 0, r G A^\ < n < M, 
and 

( ' ^/ , v) + {dT.t v) + (rft/,", dv) 

(6.6) ,(M^-M")-(n^-^-n"-i) , „ , , 

= {— ^,^') + ( nl,v)~{ptv) 

= w G A^ 1 < n < M, 

where 2;" G L^A'^ is defined by the last equation. We easily see that 

1 /•*" At 



By Theorem 13. 2[ the last term on the right hand side is bounded by Ch^\\dd*u\\Lac(^jjr-2), 
and, by (15.40 . the first term on the right hand side by 

~^\rh\\Lmt"-\t^],H^)- 

Thus we have proved: 
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Proposition 6.1. 

< ^\\Utt\\L^{L^) + Ch''{\\dd*u\\L^(^Hr-2) + -^\\Ut\\L^(lt^-i,t%Hr))- 

We shall now use equations ( 16 .Sp and (16. 6p to derive bounds on Hh and Uh in terms of z. 
Toward this end we adopt the notations 

M M 
n=l n=l 

Proposition 6.2. Assume uq G D{L). Then 

\\Uh\\lo^(L^) + \\^h\\p(L^) + \\dUh\\l2(L^) < C{\\UhiO)\\ + |k||/i(L2)), 

||Sh||/oo(i2) + ||rfS/j||/2(i2) < C(||c/;;t/ft(0)|| + \\z\\i2(^l2)). 

Proof. Take r = S^J G A{;~^ in ([63]), i; = f/^ G in ([61]), and add to obtain 

(6.7) At(||S^f + \\dUa') + r."f = (f/r^ + At^", f/,"), 
which implies 

\\un<\\u'r'\\+^t\\z-\\. 

By iteration, 

(6.8) \mi^^L2)<\K\\+Atj2m 

n=l 

which is the desired bound on Uh- To get the bound on T^h and dUh, we derive from ( 16. 7p 
that 

l\\ui:r - Imr'f + m^r + \\du;:r) < (mz-, u-) < At|i^"iiiit/,ii,.(^2). 

Summing then gives 

\\\Ulfr - \\\U'nr + l|S.||f2(,2) + \\dUh\W,,^ < ||f/.||z^(L^)|||klbi(L^), 
and so, by (16. Sp . 

3 
T 

which completes the proof of the first inequality. 

To prove the second inequality, we take r = G A^~^ in ( 16. 5 p at both time level n — 1 
and level n. This gives 

(6.9) (Sr\E5:) = (rfS«,f/r^), {K,^l) = idElUj:), l<n<M. 
Next take v = rfS^ G in ( 16. 6 p and substitute (16. 9p to get 

(S^ - Sr\S^) + AtlldS^f = A^(^^dSD, 1 < n < M, 

whence 

l|S5:iP-||Sr^f + At|ME5:|p<At||z"f. 
Again we get a telescoping sum, so 

ll^/i||/°°(L2) 

+ IMS,||f2(^2)<C(||E°|P+ 11.11^2(^2)). 
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This implies the second inequahty and so completes the proof of the proposition. □ 

As in Section [5l we choose the initial data to equal the elliptic projection of uq for 
simplicity. Then Uh{0) = and the right-hand sides of the inequalities in Proposition 16.21 
simplify. Bounding them via Proposition 16.11 we get for the first 

(6.10) ||f//i||;oo(i2) + ||S/i||i2(i2) + \\dUh\\l2(L'2) 

< CAt\\Utt\\L°-(L2) + Ch''{\\dd*u\\Loo(^H--2) + ||Mt||^i(^r)). 

For the second, we bound the t"]) norm in Proposition 16.11 by At times the L°° 



norm, and substitute the resulting bound for z in the second estimate of Proposition 16. 2[ 
obtaining 

(6.11) ||S,,||i^(i2) + ||rfS/j||,2(i2) < CAt\\uu\\L'^{L2) + Ch''{\\dd*u\\L'^(^H--2) + \\Ut\\L^(^Hr))- 

Combining fl6.10p . (16. lip with the estimates in Theorem 13.21 for the elliptic projection, we 
obtain the main result of the section. 

Theorem 6.3. Under the same assumptions as Theorem \5.3[ Let {a,u) be the solution of 
( [j.^[ ), U.5\) satisfying U.3\) . and {a^, u"^) the solution of l\6.1\) . I\6.^) with equal to the 
elliptic projection of uq. Then, we have the following error estimates for the fully discrete 
finite element method: 

\\<y - <yh\\p{L^) < CAt\\utt\\L^{L^) + Ch''{\\Ut\\Li{H-) + \\d*u\\L^{Hr)), 
\W - 0'?i||/o=(L2) < CAt\\Utt\\L°°{L2) + Ch''{\\Ut\\Loc(^Hr) + \\d*u\\Loo(^Hr)), 
\\d{cr - 0-h,)||;2(i2) < CAt\\Utt\\L^{L2) + Ch^'dlutllLoof^Hr) + \\dd*u\\Lc^(^H-)), 

\\U - n/i||,oo(i2) < CAt\\Utt\\L^{L^) + C/l''(||u||ioo(jJr) + \\Ut\\Ll(^Hr))), 

\\d{u - Uh)\\i2(^L2) < CAt\\utt\\L^(L2) + Ch''{\\ut\\Li{H-) + I Mm 1 1 (//'■) + \\dd*u\\L^(^Hr-i))- 

The error estimates are analogous to those of Theorem 15.31 for the semidiscrete solution, 
with each containing an additional 0{At) term coming from the time discretization. For 
each quantity, the error is of order 0{At + h^). 

7. Numerical examples 

In this section, we present the results of simple numerical computations verifying the 
theory above. 

First we compute a two-dimensional example for the 1-form Hodge heat equation. Using 
vector proxies, we may write the parabolic equations (ll.ll) - (ll.3p as 

Ut + (curl rot -V div)^ = / in x [0, T] , 

u ■ n = rot n = on dVt x [0, T], n( ■ , 0) = in fi, 

where 

du2 dui , , du du ^ 
rot u = — — , curl u 



dxi dx2 dx2 dxi 

We choose to be a square annulus [0, 1] x [0, 1]\[0.25, 0.75] x [0.25, 0.75] and take the exact 
solution as 

( 100x(x - l)(x - 0.25)(x - 0.75)A 
I 1002/(2/ -l)(y-0.25)(y- 0.75)tj ' 



u 
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Note that this function is not orthogonal to 1-harmonic forms on f2. We use the finite element 
spaces VrA^iTh) (Lagrange elements of degree r) for a = — divw and V~A^{Th) (Raviart- 
Thomas elements) for u, starting with an initial unstructured mesh, and then refining it 
uniformly. We take At = 0.0001 and compute the error at time T = 0.01 (after 100 time 
steps). Tables [1] and |2] show the results for r = 1 and 2 respectively. The rates of convergence 
are just as predicted by the theory. 



mesh size 


\W - <^h\\ 


rate 


l|V(a-fr/,)|| 


rate 


\\U - Uh\\ 


rate 


h 


0.0008490 


1.99 


0.1026276 


1.01 


0.0010586 


0.96 


h/2 


0.0002132 


1.99 


0.0512846 


1.00 


0.0005341 


0.99 


h/A 


0.0000534 


2.00 


0.0256528 


1.00 


0.0002678 


1.00 


h/8 


0.0000133 


2.00 


0.0128295 


1.00 


0.0001340 


1.00 


Table 1 


. Computation with PiA° x in 


two dimensions. 


mesh size 


Ik - (^hW 


rate 


l|V(a-a,)l| 


rate 


\\U - Uh\\ 


rate 


h 


0.0000093 


3.03 


0.0016510 


2.03 


0.0000705 


1.97 


h/2 


0.0000012 


3.00 


0.0004119 


2.00 


0.0000178 


1.99 


h/A 


0.0000001 


3.00 


0.0001031 


2.00 


0.0000045 


1.99 


h/8 


0.0000000 


3.04 


0.0000258 


2.00 


0.0000011 


2.00 



Table 2. Computation with P2A^ x A^ in two dimensions. 



For the second example, we let Q be the unit cube [0, 1] x [0, 1] x [0, 1] in M^, and again solve 
the 1-form Hodge heat equation. Using vector proxies, the initial-boundary value problem 
becomes 

Ut + (curl curl —V div)^ = / in x [0, T] 
u ■ n = 0, curlw x n = on dfl x [0, T], u{- ,0) = uq in Q. 
We take the exact solution to be 

(sin(7rxi)t\ 
sin(7rx2)t I . 
sin(7rx3)ty 

Table [3] shows the errors and rates of convergence for linear elements on a sequence of uniform 
meshes, again at time T = 0.01 after 100 time steps. Once again, the rates of convergence 
are just as predicted by the theory. 
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